home *** CD-ROM | disk | FTP | other *** search
/ Graphics Plus / Graphics Plus.iso / msdos / raytrace / pov / source / spheres.c < prev    next >
C/C++ Source or Header  |  1993-07-29  |  12KB  |  417 lines
  1. /****************************************************************************
  2. *                spheres.c
  3. *
  4. *  This module implements the sphere primitive.
  5. *
  6. *  from Persistence of Vision Raytracer
  7. *  Copyright 1993 Persistence of Vision Team
  8. *---------------------------------------------------------------------------
  9. *  NOTICE: This source code file is provided so that users may experiment
  10. *  with enhancements to POV-Ray and to port the software to platforms other 
  11. *  than those supported by the POV-Ray Team.  There are strict rules under
  12. *  which you are permitted to use this file.  The rules are in the file
  13. *  named POVLEGAL.DOC which should be distributed with this file. If 
  14. *  POVLEGAL.DOC is not available or for more info please contact the POV-Ray
  15. *  Team Coordinator by leaving a message in CompuServe's Graphics Developer's
  16. *  Forum.  The latest version of POV-Ray may be found there as well.
  17. *
  18. * This program is based on the popular DKB raytracer version 2.12.
  19. * DKBTrace was originally written by David K. Buck.
  20. * DKBTrace Ver 2.0-2.12 were written by David K. Buck & Aaron A. Collins.
  21. *
  22. *****************************************************************************/
  23.  
  24. #include "frame.h"
  25. #include "vector.h"
  26. #include "povproto.h"
  27.  
  28. #ifndef Sphere_Tolerance
  29. #define Sphere_Tolerance 1.0e-8
  30. #endif
  31.  
  32. METHODS Sphere_Methods =
  33.   { 
  34.   All_Sphere_Intersections,
  35.   Inside_Sphere, Sphere_Normal,
  36.   Copy_Sphere,
  37.   Translate_Sphere, Rotate_Sphere,
  38.   Scale_Sphere, Transform_Sphere, Invert_Sphere,
  39.   Destroy_Sphere
  40. };
  41.  
  42. METHODS Ellipsoid_Methods =
  43.   { 
  44.   All_Ellipsoid_Intersections,
  45.   Inside_Ellipsoid, Ellipsoid_Normal,
  46.   Copy_Sphere,
  47.   Translate_Sphere, Rotate_Sphere,
  48.   Scale_Sphere, Transform_Sphere, Invert_Sphere,
  49.   Destroy_Sphere
  50. };
  51.  
  52. extern RAY *CM_Ray;
  53. extern long Ray_Sphere_Tests, Ray_Sphere_Tests_Succeeded;
  54.  
  55. int All_Sphere_Intersections (Object, Ray, Depth_Stack)
  56. OBJECT *Object;
  57. RAY *Ray;
  58. ISTACK *Depth_Stack;
  59.   {
  60.   DBL Depth1, Depth2;
  61.   VECTOR IPoint;
  62.   register int Intersection_Found;
  63.  
  64.   Intersection_Found = FALSE;
  65.  
  66.   if (Intersect_Sphere (Ray, (SPHERE*) Object, &Depth1, &Depth2))
  67.     {
  68.     VScale (IPoint, Ray->Direction, Depth1);
  69.     VAddEq (IPoint, Ray->Initial);
  70.  
  71.     if (Point_In_Clip (&IPoint, Object->Clip))
  72.       {
  73.       push_entry(Depth1,IPoint,Object,Depth_Stack);
  74.       Intersection_Found = TRUE;
  75.       }
  76.  
  77.     if (Depth2 != Depth1)
  78.       {
  79.       VScale (IPoint, Ray->Direction, Depth2);
  80.       VAddEq (IPoint, Ray->Initial);
  81.  
  82.       if (Point_In_Clip (&IPoint, Object->Clip))
  83.         {
  84.         push_entry(Depth2,IPoint,Object,Depth_Stack);
  85.         Intersection_Found = TRUE;
  86.         }
  87.       }
  88.     }
  89.   return (Intersection_Found);
  90.   }
  91.  
  92. int All_Ellipsoid_Intersections (Object, Ray, Depth_Stack)
  93. OBJECT *Object;
  94. RAY *Ray;
  95. ISTACK *Depth_Stack;
  96.   {
  97.   DBL Depth1, Depth2, len;
  98.   VECTOR IPoint, dv;
  99.   register int Intersection_Found;
  100.   RAY New_Ray;
  101.  
  102.   /* Transform the ray into the sphere's space */
  103.   MInvTransPoint(&New_Ray.Initial, &Ray->Initial, ((SPHERE *)Object)->Trans);
  104.   MInvTransDirection(&New_Ray.Direction, &Ray->Direction, ((SPHERE *)Object)->Trans);
  105.  
  106.   VDot(len, New_Ray.Direction, New_Ray.Direction);
  107.   if (len == 0.0)
  108.     return 0;
  109.   len = 1.0 / sqrt(len);
  110.   VScaleEq(New_Ray.Direction, len);
  111.  
  112.   Intersection_Found = FALSE;
  113.  
  114.   if (Intersect_Sphere (&New_Ray, (SPHERE*) Object, &Depth1, &Depth2))
  115.     {
  116.     VScale (IPoint, New_Ray.Direction, Depth1);
  117.     VAddEq (IPoint, New_Ray.Initial);
  118.     if (((SPHERE *)Object)->Trans != NULL)
  119.       MTransPoint(&IPoint, &IPoint, ((SPHERE *)Object)->Trans);
  120.  
  121.     VSub(dv, IPoint, Ray->Initial);
  122.     VLength(len, dv);
  123.  
  124.     if (Point_In_Clip (&IPoint, Object->Clip))
  125.       {
  126.       push_entry(len,IPoint,Object,Depth_Stack);
  127.       Intersection_Found = TRUE;
  128.       }
  129.  
  130.     if (Depth2 != Depth1)
  131.       {
  132.       VScale (IPoint, New_Ray.Direction, Depth2);
  133.       VAddEq (IPoint, New_Ray.Initial);
  134.       if (((SPHERE *)Object)->Trans != NULL)
  135.         MTransPoint(&IPoint, &IPoint, ((SPHERE *)Object)->Trans);
  136.  
  137.       VSub(dv, IPoint, Ray->Initial);
  138.       VLength(len, dv);
  139.  
  140.       if (Point_In_Clip (&IPoint, Object->Clip))
  141.         {
  142.         push_entry(len,IPoint,Object,Depth_Stack);
  143.         Intersection_Found = TRUE;
  144.         }
  145.       }
  146.     }
  147.   return (Intersection_Found);
  148.   }
  149.  
  150. int Intersect_Sphere (Ray, Sphere, Depth1, Depth2)
  151. RAY *Ray;
  152. SPHERE *Sphere;
  153. DBL *Depth1, *Depth2;
  154.   {
  155.   VECTOR Origin_To_Center;
  156.   DBL OCSquared, t_Closest_Approach, Half_Chord, t_Half_Chord_Squared;
  157.   short inside;
  158.  
  159.   Ray_Sphere_Tests++;
  160.   if (Ray == CM_Ray) 
  161.     {
  162.     if (!Sphere->CMCached) 
  163.       {
  164.       VSub (Sphere->CMOtoC, Sphere->Center, Ray->Initial);
  165.       VDot (Sphere->CMOCSquared, Sphere->CMOtoC, Sphere->CMOtoC);
  166.       Sphere->CMinside = (Sphere->CMOCSquared < Sphere->Radius_Squared);
  167.       Sphere->CMCached = TRUE;
  168.       }
  169.     VDot (t_Closest_Approach, Sphere->CMOtoC, Ray->Direction);
  170.     if (!Sphere->CMinside && (t_Closest_Approach < Sphere_Tolerance))
  171.       return (FALSE);      
  172.     t_Half_Chord_Squared = Sphere->Radius_Squared - Sphere->CMOCSquared +
  173.     (t_Closest_Approach * t_Closest_Approach);
  174.     }
  175.   else 
  176.     {
  177.     VSub (Origin_To_Center, Sphere->Center, Ray->Initial);
  178.     VDot (OCSquared, Origin_To_Center, Origin_To_Center);
  179.     inside = (OCSquared < Sphere->Radius_Squared);
  180.     VDot (t_Closest_Approach, Origin_To_Center, Ray->Direction);
  181.     if (!inside && (t_Closest_Approach < Sphere_Tolerance))
  182.       return (FALSE);
  183.  
  184.     t_Half_Chord_Squared = Sphere->Radius_Squared - OCSquared +
  185.     (t_Closest_Approach * t_Closest_Approach);
  186.     }
  187.  
  188.   if (t_Half_Chord_Squared < Sphere_Tolerance)
  189.     return (FALSE);
  190.  
  191.   Half_Chord = sqrt (t_Half_Chord_Squared);
  192.   *Depth1 = t_Closest_Approach + Half_Chord;
  193.   *Depth2 = t_Closest_Approach - Half_Chord;
  194.  
  195.   if ((*Depth1 < Sphere_Tolerance) || (*Depth1 > Max_Distance))
  196.     if ((*Depth2 < Sphere_Tolerance) || (*Depth2 > Max_Distance))
  197.       return (FALSE);
  198.     else
  199.       *Depth1 = *Depth2;
  200.   else
  201.     if ((*Depth2 < Sphere_Tolerance) || (*Depth2 > Max_Distance))
  202.       *Depth2 = *Depth1;
  203.  
  204.   Ray_Sphere_Tests_Succeeded++;
  205.   return (TRUE);
  206.   }
  207.  
  208. int Inside_Sphere (IPoint, Object)
  209. VECTOR *IPoint;
  210. OBJECT *Object;
  211.   {
  212.   VECTOR Origin_To_Center;
  213.   DBL OCSquared;
  214.  
  215.   VSub (Origin_To_Center, ((SPHERE *)Object)->Center, *IPoint);
  216.   VDot (OCSquared, Origin_To_Center, Origin_To_Center);
  217.  
  218.   if (((SPHERE *)Object)->Inverted)
  219.     return (OCSquared - ((SPHERE *)Object)->Radius_Squared > Sphere_Tolerance);
  220.   else
  221.     return (OCSquared - ((SPHERE *)Object)->Radius_Squared < Sphere_Tolerance);
  222.   }
  223.  
  224. int Inside_Ellipsoid (IPoint, Object)
  225. VECTOR *IPoint;
  226. OBJECT *Object;
  227.   {
  228.   VECTOR Origin_To_Center;
  229.   DBL OCSquared;
  230.   VECTOR New_Point;
  231.  
  232.   /* Transform the point into the sphere's space */
  233.   MInvTransPoint(&New_Point, IPoint, ((SPHERE *)Object)->Trans);
  234.   VSub (Origin_To_Center, ((SPHERE *)Object)->Center, New_Point);
  235.   VDot (OCSquared, Origin_To_Center, Origin_To_Center);
  236.  
  237.   if (((SPHERE *)Object)->Inverted)
  238.     return (OCSquared - ((SPHERE *)Object)->Radius_Squared > Sphere_Tolerance);
  239.   else
  240.     return (OCSquared - ((SPHERE *)Object)->Radius_Squared < Sphere_Tolerance);
  241.   }
  242.  
  243. void Sphere_Normal (Result, Object, IPoint)
  244. OBJECT *Object;
  245. VECTOR *Result, *IPoint;
  246.   {
  247.   VSub (*Result, *IPoint, ((SPHERE *)Object)->Center);
  248.   VScaleEq (*Result, ((SPHERE *)Object)->Inverse_Radius);
  249.   }
  250.  
  251. void Ellipsoid_Normal (Result, Object, IPoint)
  252. OBJECT *Object;
  253. VECTOR *Result, *IPoint;
  254.   {
  255.   VECTOR New_Point;
  256.  
  257.   /* Transform the point into the sphere's space */
  258.   MInvTransPoint(&New_Point, IPoint, ((SPHERE *)Object)->Trans);
  259.  
  260.   VSub (*Result, New_Point, ((SPHERE *)Object)->Center);
  261.   VScaleEq (*Result, ((SPHERE *)Object)->Inverse_Radius);
  262.  
  263.   MTransNormal(Result, Result, ((SPHERE *)Object)->Trans);
  264.   VNormalize(*Result, *Result);
  265.   }
  266.  
  267. void *Copy_Sphere (Object)
  268. OBJECT *Object;
  269.   {
  270.   SPHERE *New;
  271.  
  272.   New = Create_Sphere ();
  273.   *New = *((SPHERE *) Object);
  274.  
  275.   New->Trans = Copy_Transform(((SPHERE *)Object)->Trans);
  276.   return (New);
  277.   }
  278.  
  279. void Translate_Sphere (Object, Vector)
  280. OBJECT *Object;
  281. VECTOR *Vector;
  282.   {
  283.   TRANSFORM Trans;
  284.  
  285.   if (((SPHERE *)Object)->Trans == NULL)
  286.     {
  287.     VAddEq (((SPHERE *) Object)->Center, *Vector);
  288.     VAddEq(Object->Bounds.Lower_Left, *Vector);
  289.     }
  290.   else
  291.     {
  292.     Compute_Translation_Transform(&Trans, Vector);
  293.     Transform_Sphere(Object, &Trans);
  294.     }
  295.   }
  296.  
  297. void Rotate_Sphere (Object, Vector)
  298. OBJECT *Object;
  299. VECTOR *Vector;
  300.   {
  301.   TRANSFORM Trans;
  302.   SPHERE *Sphere = (SPHERE *) Object;
  303.  
  304.   Compute_Rotation_Transform (&Trans, Vector);
  305.  
  306.   if (Sphere->Trans == NULL)
  307.     {
  308.     MTransPoint(&Sphere->Center, &Sphere->Center, &Trans);
  309.     Make_Vector(&Sphere->Bounds.Lower_Left,
  310.       Sphere->Center.x - Sphere->Radius,
  311.       Sphere->Center.y - Sphere->Radius,
  312.       Sphere->Center.z - Sphere->Radius);
  313.     Make_Vector(&Sphere->Bounds.Lengths,
  314.       2.0 * Sphere->Radius,
  315.       2.0 * Sphere->Radius,
  316.       2.0 * Sphere->Radius);
  317.     }
  318.   else
  319.     {
  320.     Transform_Sphere (Object, &Trans);
  321.     }
  322.   }
  323.  
  324. void Scale_Sphere (Object, Vector)
  325. OBJECT *Object;
  326. VECTOR *Vector;
  327.   {
  328.   SPHERE *Sphere = (SPHERE *) Object;
  329.   TRANSFORM Trans;
  330.  
  331.   if ((Vector->x != Vector->y) || (Vector->x != Vector->z))
  332.     if (Sphere->Trans == NULL)
  333.       {
  334.       Sphere->Methods = &Ellipsoid_Methods;
  335.       Sphere->Trans = Create_Transform();
  336.       }
  337.  
  338.   if (Sphere->Trans == NULL)
  339.     {
  340.     VScaleEq (Sphere->Center, Vector->x);
  341.     Sphere->Radius *= Vector->x;
  342.     Sphere->Radius_Squared = Sphere->Radius * Sphere->Radius;
  343.     Sphere->Inverse_Radius = 1.0 / Sphere->Radius;
  344.     Make_Vector(&Sphere->Bounds.Lower_Left,
  345.       Sphere->Center.x - Sphere->Radius,
  346.       Sphere->Center.y - Sphere->Radius,
  347.       Sphere->Center.z - Sphere->Radius);
  348.     Make_Vector(&Sphere->Bounds.Lengths,
  349.       2.0 * Sphere->Radius,
  350.       2.0 * Sphere->Radius,
  351.       2.0 * Sphere->Radius);
  352.     }
  353.   else
  354.     {
  355.     Compute_Scaling_Transform(&Trans, Vector);
  356.     Transform_Sphere(Object, &Trans);
  357.     }
  358.   }
  359.  
  360. void Invert_Sphere (Object)
  361. OBJECT *Object;
  362.   {
  363.   ((SPHERE *) Object)->Inverted ^= TRUE;
  364.   }
  365.  
  366. SPHERE *Create_Sphere ()
  367.   {
  368.   SPHERE *New;
  369.  
  370.   if ((New = (SPHERE *) malloc (sizeof (SPHERE))) == NULL)
  371.     MAError ("sphere");
  372.  
  373.   INIT_OBJECT_FIELDS(New, SPHERE_OBJECT, &Sphere_Methods)
  374.     Make_Vector (&(New->Center), 0.0, 0.0, 0.0);
  375.   New->Radius = 1.0;
  376.   New->Radius_Squared = 1.0;
  377.   New->Inverse_Radius = 1.0;
  378.   New->CMCached = FALSE;
  379.   New->Trans = NULL;
  380.   New->Inverted = FALSE;
  381.   /* CMOtoC, CMOtoCSquared and CMinside are only valid when CMCached */
  382.   return (New);
  383.   }
  384.  
  385. void Transform_Sphere (Object, Trans)
  386. OBJECT *Object;
  387. TRANSFORM *Trans;
  388.   {
  389.   SPHERE *Sphere = (SPHERE *)Object;
  390.  
  391.   if (Sphere->Trans == NULL)
  392.     {
  393.     Sphere->Methods = &Ellipsoid_Methods;
  394.     Sphere->Trans = Create_Transform();
  395.     }
  396.  
  397.   Compose_Transforms(Sphere->Trans, Trans);
  398.  
  399.   Make_Vector(&Sphere->Bounds.Lower_Left,
  400.     Sphere->Center.x - Sphere->Radius,
  401.     Sphere->Center.y - Sphere->Radius,
  402.     Sphere->Center.z - Sphere->Radius);
  403.   Make_Vector(&Sphere->Bounds.Lengths,
  404.     2.0 * Sphere->Radius,
  405.     2.0 * Sphere->Radius,
  406.     2.0 * Sphere->Radius);
  407.  
  408.   recompute_bbox(&Object->Bounds, Sphere->Trans);
  409.   }
  410.  
  411. void Destroy_Sphere (Object)
  412. OBJECT *Object;
  413.   {
  414.   Destroy_Transform(((SPHERE *)Object)->Trans);
  415.   free (Object);
  416.   }
  417.